Geosci. Model Dev., 5, 683-707, 2012 
www.geosci-model-dev.net/ 5/683/2012/ 
doi:10.5194/gmd-5-683-2012 
© Author(s) 2012. CC Attribution 3.0 License. 




Geoscientific 
Model Development 


The role of phytoplankton dynamics in the seasonal and interannual 
variability of carbon in the subpolar North Atlantic - a modeling 
study 

S. R. Signorini 1 ’ 9 , S. Hakkinen 2 , K. Gudmundsson 3 , A. Olsen 4 ’ 5 *, A. M. Omar 4 5 , J. Olafsson 6 , G. Reverdin 7 , 

S. A. Henson 8 , C. R. McClain 9 , and D. L. Worthen 210 

'Science Applications International Corp., Beltsville, Maryland, USA 

2 NASA Goddard Space Flight Center, Cryospheric Sciences Branch, Greenbelt, Maryland, USA 

3 Marine Research Institute, RO. Box 1390, 121, Reykjavik, Iceland 

4 Uni Bjerknes Centre, Allegaten 55, 5007 Bergen, Norway 

5 Bjerknes Centre for Climate Research, Allegaten 55, 5007 Bergen, Norway 

6 University of Iceland and Marine Research Institute, Reykjavik, Iceland 

7 Laboratoire d’Oceanographie Dynamique et de Climatologie, IPSL, Boite 100, 4, Place Jussieu 75252, Paris 
8 National Oceanography Centre, European Way, Southampton, SOM 3ZH, UK 
9 NASA Goddard Space Flight Center, Ocean Ecology Branch, Greenbelt, Code 614.2, Maryland, USA 
*°Wyle Information Systems Group, McLean, Virginia, USA 

*now at: Institute of Marine Research, P.O. Box 1870 Nordnes, 5817 Bergen, Norway 
Correspondence to: S. R. Signorini (sergio.signorini@nasa.gov) 

Received: 3 January 2011 - Published in Geosci. Model Dev. Discuss.: 3 February 2011 
Revised: 22 March 2012 - Accepted: 29 March 2012 - Published: 15 May 2012 


Abstract. We developed an ecosystem/biogeochemical 
model system, which includes multiple phytoplankton func- 
tional groups and carbon cycle dynamics, and applied it to in- 
vestigate physical-biological interactions in Icelandic waters. 
Satellite and in situ data were used to evaluate the model. 
Surface seasonal cycle amplitudes and biases of key parame- 
ters (DIC, TA, pC0 2 , air-sea C0 2 flux, and nutrients) are sig- 
nificantly improved when compared to surface observations 
by prescribing deep water values and trends, based on avail- 
able data. The seasonality of the coccolithophore and “other 
phytoplankton” (diatoms and dinoflagellates) blooms is in 
general agreement with satellite ocean color products. Nutri- 
ent supply, biomass and calcite concentrations are modulated 
by light and mixed layer depth seasonal cycles. Diatoms 
are the most abundant phytoplankton, with a large bloom 
in early spring and a secondary bloom in fall. The diatom 
bloom is followed by blooms of dinoflagellates and coccol- 
ithophores. The effect of biological changes on the seasonal 
variability of the surface ocean pCOi is nearly twice the tem- 
perature effect, in agreement with previous studies. The in- 
clusion of multiple phytoplankton functional groups in the 


model played a major role in the accurate representation of 
C0 2 uptake by biology. For instance, at the peak of the 
bloom, the exclusion of coccolithophores causes an increase 
in alkalinity of up to 4 nmol kg 1 with a corresponding in- 
crease in DIC of up to 16 pmolkg"'. During the peak of the 
bloom in summer, the net effect of the absence of the coc- 
colithophores bloom is an increase in pCOj of more than 
20 patm and a reduction of atmospheric C0 2 uptake of more 
than 6mmolm _2 d _1 . On average, the impact of coccol- 
ithophores is an increase of air-sea C0 2 flux of about 27 %. 
Considering the areal extent of the bloom from satellite im- 
ages within the Irminger and Icelandic Basins, this reduction 
translates into an annual mean of nearly 1500 tonnes Cyr - '. 


1 Introduction 

Coupled biogeochemical-physical numerical models to- 

gether with observations are an essential tool to understand 

the interaction between physical and biological processes 

that create the observed temporal variability, on time scales 
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ranging from days to years. Here, we developed a coupled 
biogeochemical-physical model system, which includes mul- 
tiple phytoplankton functional groups and carbon cycle dy- 
namics, and applied it to assess the role of phytoplankton 
dynamics in the seasonal and interannual variability of car- 
bon in the subpolar North Atlantic. The rationale for the site 
selection relates to a combination of environmental and bio- 
logical factors, e.g., large range in mixed layer depth (MLD), 
low mean horizontal advection, clear seasonal succession of 
phytoplankton species, and bloom intensity. Spring phy- 
toplankton blooms in the surface water of the North At- 
lantic Ocean and adjacent seas are known to cause a pre- 
cipitous reduction of surface water partial pressure of car- 
bon dioxide (pC02) and the concentrations of CO 2 and nu- 
trients on a seasonal basis (Takahashi et al., 1993). In ad- 
dition, research during the last decade has increased aware- 
ness of the relationship between key phytoplankton groups 
and their pivotal roles in the biogeochemical cycles of a 
range of elements (Boyd et ah, 2010). For example, coccol- 
ithophores occupy a central role in the carbon cycle by the 
conversion of dissolved inorganic carbon (DIC) to both par- 
ticulate organic carbon (POC) and particulate inorganic car- 
bon (PIC) forms, albeit uncertainties exist over whether coc- 
colithophore blooms are net sinks or sources of CO 2 to the at- 
mosphere (Boyd and Trull, 2007). The ratio of PIC to POC in 
exported biogenic matter determines the relative strength of 
the biological carbon pump and consequently the flux of CO 2 
across the surface ocean-atmosphere interface. Other func- 
tional groups, such as diatoms and dinoflagellates, are also 
important in the sequestration of carbon from the atmosphere 
to the deep ocean. The ability of some species of diatoms 
to form chains with built-in silica ballast and to produce 
large fast-sinking aggregates during the declining phase of 
the blooms means that they are key vectors in exporting and 
sequestering POC to the deep ocean (Lampitt, 1985). The 
competitive success of phytoplankton functional groups has 
important links to the environmental conditions. Although 
the predictive outcome of competition between two or more 
species has seldom been tested in the marine environment, 
some laboratory and mesocosm experiments reveal interest- 
ing results. Elaborate laboratory culture competition exper- 
iments with multiple marine phytoplankton species (coccol- 
ithophores excluded) (Sommer, 1994) and mesocosm experi- 
ments (Egge and Aksnes, 1992) show that diatoms dominate 
at high Si:N ratios, while diverse flagellates tend to dominate 
at low Si:N ratios. As Si is not regenerated above the MLD 
during a season, the growth period for diatoms fades out as 
the pycnoline is enforced and may be revoked every time a 
vertical mixing event occurs, as long as the available light 
is sufficient. Sommer (1994) showed no apparent effect of 
light intensities on the succession of functional groups, al- 
though it may exert some influence at the species level on 
the transition along the Si:N gradient. Nevertheless, light in- 
tensity has been shown to be very important on the blooms 
of certain species of coccolithophores. The large E. hux- 
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leyi bloom south of Iceland in 1991 was modeled and high 
light was shown to be a positive cause (Tyrrell and Taylor, 
1996). Under the intense bloom conditions , mixing was not 
very deep and therefore average light intensities in the sur- 
face layer were high. 

The physical characteristics of the North Atlantic subpolar 
waters also have an impact on the biogeochemical and carbon 
cycle variability. The waters surrounding Iceland are charac- 
terized by the cold polar water of the East Greenland Current 
and Arctic water of the East Icelandic Current from the north, 
and the warm North Atlantic water of the Irminger Current 
from the south (Gudmundsson, 1998). Figure la shows a 
map of the subpolar North Atlantic and Nordic Seas and the 
mean location of the Arctic Front (AF). Sea ice is present 
throughout the year along the east coast of Greenland, with 
the ice edge extending farthest offshore in winter-spring, and 
retreating in the fall. Ice melt in summer-fall freshens (32- 
34psu) the surface waters north of the AF. South of the AF 
waters are saltier (~35psu) and warmer (8 to 12 °C), as a 
result of Atlantic water intrusions. These hydrographic char- 
acteristics have an impact on the seasonal vertical mixing, 
as shown in model simulations. During winter and spring, 
MLDs average 100 to 500 m south of the AF, with deeper 
values in the Irminger and Icelandic Basins; in summer, the 
vertical stratification is significant, with MLDs less than 20 m 
(de Boyer Montegut et al., 2004; de Boisseson et al., 2010; 
Carton et al., 2008). North of the AF, the MLDs are shallow 
(< 40 m) throughout the year due to the southward advec- 
tion of fresher polar and Arctic waters, and ice melting dur- 
ing summer-fall. Some accounts of winter deep mixing are 
available for the region based on regular XBT surveys and 
CTD casts (G. Reverdin, personal communication, 2012). 
MLDs as deep as 500 m were found in the Irminger Current, 
in particular west of the Reykjanes Ridge, and occasionally 
close to Greenland and southwest of Iceland. In the central 
Irminger Sea, MLDs are not as deep, mostly ranging from 
100 to 150 m. However, there are indications from Argos 
floats that deep convective chimneys occasionally occur in 
the central Irminger Sea (as in early 2008) down to nearly 
1000 m (Vage et al., 2009). These large seasonal changes 
in stratification and vertical mixing play an important role in 
the euphotic zone nutrient renewal and on the onset and dura- 
tion of the phytoplankton spring bloom (Henson et al., 2006). 
Variability in the intensity of primary production in general 
and of the timing of spring bloom in particular (Henson et 
al., 2009) affects the population dynamics of higher trophic 
levels, such as the commercially important Atlanto-Scandian 
herring (Jakobsson, 1978) in the region. The annual phy- 
toplankton carbon production in the waters around Iceland is 
quite high. According to a recent estimate (Zhai et al., 2012), 
based on available regional measurements on the photosyn- 
thetic capacity and satellite biomass data, the annual primary 
production is 179 ±36 and 238±22gCm _2 y _1 for the Arc- 
tic and the Atlantic waters, respectively. In addition, sea- 
sonal and interannual changes in phytoplankton production 
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Fig. la. Map of subpolar North Atlantic and Nordic Seas. The black thick line is the Arctic Front represented by the 4 °C isotherm (from 
GPOMZ 1948-2009 climatology) and the thinner lines are the bathymetry contours (500, 1500, and 3000m). The black triangle indicates 
the site chosen for the 1-D ecosystem-carbon model. The blue circles (1994) and red crosses (2004) indicate the location of the two repeat 
CARINA transects. 



oSuratlant T. S, Nutrients, DIC, TA. pC0 2 * Mar Prod PP 
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Fig. lb. Map showing the spatial distribution of in situ PP, nutrients, 
and carbon data used for model evaluation. The rectangular box 
bounds the CPR standard area B6, from which taxonomy data were 
obtained for evaluation. 

are tightly coupled to atmospheric CO 2 uptake and surface 
ocean pCOj variability (Corbiere et ah, 2007), and therefore 
a major component of the carbon cycle. 


2 Model description 

2.1 Ice-ocean model (GPOMZ) 

Physical forcing for the 1-D ecosystem/carbon model, 
ECOlD-2.0(location shown by the black triangle in Fig. la 
at 30° W and 60° N), is provided by a three-dimensional (3- 
D) coupled ice-ocean model of the Arctic and North At- 
lantic (Mellor et ah, 2002). The sigma-coordinate version 
has been previously used to simulate North Atlantic and Arc- 
tic variability (Hakkinen, 1999, 2001; Hakkinen and Geiger 
2000). Monthly temperature and salinity profiles obtained 
from GPOMZ are interpolated to the time steps and depths 
required by ECO1D-2.0. These profiles are assimilated at 
each time step by the 1-D model. A detailed description of 
GPOMZ is provided in Appendix A. In this subsection we 
present a few GPOMZ products to describe some of the phys- 
ical properties of the region. In Appendix A we show three 
figures with physical properties derived from GPOMZ. Fig- 
ure Al shows the annual SST composites for 2009, derived 
from GPOMZ and MODIS-Aqua; Fig. A2 shows a compari- 
son between GPOMZ and MODIS SST at ECO1D-2.0 sim- 
ulation site, and Fig. A3 shows a comparison between the 
seasonal MLD from de Boyer Montegut et ah (2004) and 
GPOMZ at the ECO1D-2.0 site. Although there is a rea- 
sonable agreement between model and MODIS SST fields 
in Fig. Al, the Arctic Front is further away from the study 
site compared to MODIS data. However, as it will be shown 
later in section 4 and Appendix A, this has little or no effect 
on our results, since the bias-corrected SST at the modeling 
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Fig. 2. Seasonal maps of GPOMZ SST (left panel, a through d) 
and SSS (right panel, a through d) for 2005. The 3D ice concentra- 
tion contours (20 %, 40 %, and 60 % purple contours) derived from 
monthly NOAA/ESRL/PSD optimally interpolated data (OI) are su- 
perposed. The ECO1D-2.0 site is shown by the black triangle and 
the thick black lines indicate the Arctic Front. 


site agrees well with the data. The bias correction was done 
based on WOA05 monthly climatology and the resulting SST 
time series at the ECO1D-2.0 site was validated based on in 
situ observations and MODIS data (see caption in Fig. A2 for 
explanation). 

Figure lb shows a smaller-scale map containing the lo- 
cations of the in situ data available for model evaluation. 
The 1-D model location lies on the western edge of the 
Reykjanes Ridge and within the Irminger Current, in an 
area where winter mixing is extremely vigorous (mixed layer 
depths approaching 600 m) (Bailey et al., 2005) and spring- 
summer phytoplankton blooms are substantial (Fernandez et 
al., 1993; Flolligan et al., 1993b; Weeks et al., 1993). The 


convergence of the polar and North Atlantic water masses 
form the Arctic Front, which varies slightly in location sea- 
sonally (Fig. 2) and is identified by strong gradients of sea 
surface temperature (SST) and salinity (SSS). The AF runs 
east-west, approximately along 66° N to the west of Iceland, 
and north-south between 12° W-8° W, and then further in 
northerly direction from approximately 69° N, northeast of 
Iceland. The above features are clearly shown in Fig. 2, 
which shows the seasonal SST and SSS fields for 2005 ob- 
tained from GPOMZ (Mellor et al., 2002). GPOMZ products 
are available at monthly intervals for 1948-2009. We chose 
to show fields for 2005 as an example because it coincides 
with the year during which most of the surface ocean pCOi 
data are available (Chierici et al., 2009; Olsen et al., 2008). 

The climatologic (1954-2008) seasonal surface currents 
from GPOMZ are shown in Fig. 3. Although there may 
be some differences between the model currents and what 
is known of the regional circulation patterns (Poulain et al., 
1996; Gudmundsson, 1998), Fig. 3 shows a general cyclonic 
circulation in the Irminger Sea north of 59° N, and a conver- 
gence south of the Denmark Strait between the cyclonic cir- 
culation and the East Greenland Current, in agreement with 
observations. 

2.2 Ecosystem-carbon model (ECO1D-2.0) 

The mixed layer component of ECO1D-2.0 originates from 
an existing one-dimensional physical-biogeochemical model 
(Signorini et al., 2001a; Signorini et al., 2001b), which uti- 
lizes a turbulence closure mixed layer scheme (Mellor and 
Yamada, 1982), identical to the scheme used for GPOMZ. 
ECO1D-2.0 has a vertical coordinate system that provides 
parameter values, including horizontal velocity components, 
temperature, salinity, and the vertical diffusivity coefficient, 
k v , at each time step and grid point. The surface boundary 
layer is resolved more accurately by using a stretched verti- 
cal coordinate with higher resolution near the surface. 

The biogeochemical component of ECO1D-2.0 derives 
from a previous ecosystem model, ECO1D-1.0, configured 
for and applied at the Bermuda Atlantic Time Series (BATS) 
site, described in Signorini et al. (2003). ECO1D-2.0 in- 
cludes additional conservation equations for diatoms, coccol- 
ithophores, calcium carbonate (CaCCb), silicate (SiOi) and 
alkalinity. Figure 4 shows a diagram illustrating the model 
components and their couplings. The details of the ecosys- 
tem model are described in Appendix B. 

Figure 5 (to be discussed later in Sect. 4) shows the 
satellite-derived net primary production (NPP) and calcifica- 
tion rate (PpcaCCw) for June 1998, during which the strongest 
coccolithophore bloom occurred during the SeaWiFS mis- 
sion (Raitsos et al., 2006), with the location of the model- 
ing site (black triangle at 30° W and 60° N). NPP was ob- 
tained from monthly CbPM files and the PpcaC03 monthly 
composite was computed using the algorithm of Balch et 
al. (2007). The site lies at the western edge of the Reykjanes 
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Fig. 3. Map showing the 3-D model seasonal climatologic (1955-2008) surface currents and the location of the 1-D ecosystem model 
simulations (black triangle). 



Fig. 4. Diagram of main components of the one-dimensional bio- 
geochemical model. The color code for some of the variables in- 
dicates unique association of processes, e.g., ammonium nitrifica- 
tion, uptake of silica by diatoms only, and calcification by coccol- 
ithophores. 

Ridge in the Irminger Sea, a region of high phytoplankton 
productivity. The local depth is ~ 1500 m, but the model 
vertical grid extends to 1000 m only, which accommodates 


the deepest MLD. ECO1D-2.0 can be forced by winds and 
other relevant atmospheric parameters to calculate heat and 
freshwater fluxes originating from NCEP-II Reanalysis prod- 
ucts, which determined the 28-yr period of simulation ( 1981— 
2008). As an option, surface temperature and salinity can be 
specified at the surface instead of heat and freshwater fluxes. 
The physical component of ECO1D-2.0 is configured on a 
vertically-stretched logarithmic grid, while the biogeochem- 
ical component uses a uniform vertical grid with 1 m reso- 
lution. Deep water (z > 200 m) nutrient values are nudged 
within the bounds of the observed monthly climatology from 
the World Ocean Atlas 2005 (see details in Appendix B), and 
temperature and salinity values at all depths are nudged to 
values derived from the GPOMZ using the Newtonian re- 
laxation method (see Appendix B for details). The Newto- 
nian relaxation method (nudging) is a simple form of data 
assimilation. Surface salinity (SSS) from GPOMZ is im- 
posed at the top layer of the 1-D model, which accounts for 
all processes that alter SSS, including ice melting-freezing 
and precipitation-evaporation. The Reynolds and Smith op- 
timally interpolated (RSOI) SST is imposed at the model sur- 
face instead of heat flux. 

Deep water dissolved inorganic carbon (DIC) is nudged 
(relaxation time of 2 days) to values obtained from a data- 
derived equation (DIC vs. T) and model temperature at each 
time step (0.5 h), plus has a superimposed DIC decadal trend 
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Fig. 5. Satellite-derived net primary production (a) and calcification rate (b) for June 1998 (units in mg Cm 2 d *). The white triangle 
indicates the position of the 1-D ecosystem-carbon model simulations. Ice concentrations greater than 10% are masked in white. 
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Fig. 6. Taylor diagram showing normalized standard deviations 
(a*) and correlation coefficients (R) for model versus reference 
values, consisting of all measured variables from different available 
sources. All comparisons were done using surface values, except 
for PP which is a vertically integrated value. The black circle shows 
the result that would be obtained for a perfect fit between model and 
data, e.g., identical standard deviations and R— 1. 


consistent with observations . The total alkalinity (TA) is 
nudged in deep water to values also obtained from a data- 
derived equation (TA vs. T and S), with a relaxation time 
scale of 10 days. The DIC and TA equations are based on 
CARINA observations, with details described in the Ap- 
pendix B. 

The adopted method of forcing the one-dimensional bio- 
geochemical model (ECO1D-2.0) with physical fields from 
GPOMZ relies on the assumption that vertical mixing pro- 


cesses, e.g., convective overturning and wind mixing, play 
a much more important role in the biological variability 
than the lateral advective processes. To better quantify and 
justify this assumption, an assessment of the nitrate trans- 
port via vertical mixing and horizontal advection was per- 
formed using a combination of monthly climatologic ni- 
trate profiles from WOA05 and monthly climatologic (1981- 
2008) physical parameters (eddy vertical diffusivity, k v , 
and north and east velocity components, u and v) from 
GPOMZ. The calculation was done within ±2° of the mod- 
eling site at 60° N and 30° W, and within the depth inter- 
val from 100 to 500 m, in the upper part of the water col- 
umn. The minimum depth of 100 m was chosen to avoid 
the euphotic depth, within which nitrate is consumed by 
phytoplankton, while the maximum depth of 500 m is dic- 
tated by the availability of monthly WOA05 nitrate. Spatial 
and temporal averages were done to the resulting horizon- 
tal divergence and vertical mixing components of the nitrate 
transport, w3N03/3x + n9N03/fi(y and 3/3z (AySNC^/Sz), 
respectively. The results demonstrate that our assump- 
tion is a valid one. The annual mean local horizon- 
tal divergence is 0.005 mmol m -3 d -1 , while the vertical 
eddy diffusivity transport is 2.628 mmolm -3 d -1 , a fac- 
tor of ~550 times larger than the local horizontal diver- 
gence. The values for summer-fall (January-March and 
October-December) and spring-summer (April-September) 
are 0.007 and 0.003 mmol m -3 d -1 for horizontal diver- 
gence, and5.124 and 0.133 mmol m -3 d -1 for vertical diffu- 
sion, respectively. Although the vertical diffusion of nitrate is 
significantly reduced in spring-summer, it is still ~52 times 
larger than the horizontal divergence. 


Geosci. Model Dev., 5, 683-707, 2012 


www.geosci-model-dev.net/5/683/2012/ 


S. R. Signorini et al.: The role of phytoplankton dynamics 


689 




420 

400 

380 

360 

340 

320 

300 

280 

25 

20 

15 

10 

5 

0 

-5 

-10 




Fig. 7a. Climatologic seasonal variability of surface ocean DIC (a), TA (b), pCC> 2 , and air-sea CO2 flux for three different model runs: ref- 
erence (Experiment A, green), no DIC and TA relaxation and without bottom DIC trend (Experiment B, red), and abiotic ocean (Experiment 
C, blue). 


3 Satellite and in situ data sets 

We rely on a combination of satellite data and field obser- 
vations to provide an evaluation framework for the model. 
The satellite-derived data sets consist of SeaWiFS Chi a, 
Reynolds and Smith optimally interpolated (RSOI) SST 
(Reynolds and Smith, 1995) and sea ice concentration. These 
data were obtained from NOAA_OI_SST_V2 data provided 
by the NOAA/OAR/ESRL PSD, Boulder, Colorado, USA, 
from their web site at http://www.esrl.noaa.gov/psd/. We 
use primary production from the carbon-based productiv- 
ity model (CbPM, (Behrenfeld et al., 2005)). The satellite- 
derived PP data were obtained from the ocean productiv- 
ity web site at Oregon State University (http ://www. science. 
oregonstate.edu/ocean.productivity/). The inputs for CbPM 
are Chi a and backscatter at 443 nm , estimated using the 
Garver-Siegel-Maritorena (GSM) semi-analytical algorithm 
(Maritorena et al., 2002) and monthly SeaWiFS water leav- 
ing radiances, SeaWiFS cloud-corrected PAR data, SeaW- 
iFS mixed layer light attenuation coefficients at 490 nm, ad- 
vanced very high resolution (AVFIRR) SST, and monthly 
mean regional mixed layer depths from the Fleet Numeric 
Meteorology and Oceanography Center (FNMOC) (Mon- 
terey, California). In addition, a limited amount of in situ 
PP data was obtained from C-14 incubations conducted dur- 
ing early August 2002 on a Marine Productivity (MarProd) 
cruise sponsored by the Natural Environment Research 
Council (UK) onboard the RRS Discovery. Taxonomic data 
(cell counts) were obtained from the Continuous Plankton 


Recorder (CPR) database (standard area B6 south of Iceland, 
http://www.sahfos.ac.uk/data-archive/standard-areas.aspx). 

Surface ocean pCOi data for 2005 were acquired from 
onboard the container ship M/V Nuka Arctica (Chierici et al., 
2009; Olsen et al., 2008). In addition, we used DIC, alka- 
linity, surface ocean pCOj, and nutrient in situ data from 
SURATLANTE (Corbiere et al., 2007; Metzl et al., 2010), 
and data from the CARbon IN the Atlantic Ocean (CA- 
RINA) Data Synthesis Project (Key et al. (2010); http://www. 
earth-syst-sci-data.net/2/105/2010/essd-2- 105-2010.pdf). 

The locations of all in situ data used to validate the 1-D 
ecosystem model are shown in Figure lb. The atmospheric 
pC0 2 required to obtain the CO 2 flux at the atmosphere- 
ocean interface was obtained from GLOBALVIEW-CO 2 
(NOAA, ftp.cmdl.noaa.gov). Note that these data sets used 
for model validation are from the surface or near-surface, 
which are not part of the deep water data utilized for model 
forcing and boundary conditions. 


4 Ecosystem-carbon model evaluation 

Using the data sets shown in Fig. lb and discussed in Sect. 3, 
model versus data comparisons were compiled and analyzed. 
Due to the nature of the spatial distribution of the in situ 
data surrounding the model location, it is expected that some 
of the mismatches between model and observations are in- 
herently related to the spatial variability and patchiness of 
measured quantities. For instance, the patchiness of primary 
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Fig. 7b. Time series of surface DIC (top) and ocean pCCb (bottom) from two model experiments, one including deep water DIC and TA 
relaxation and bottom DIC trend (Experiment A, black lines) and one without it (Experiment B, purple lines). The red, blue and black circles 
represent the available data for comparison. 


production and calcite production can be readily seen in the 
satellite maps shown in Fig. 5. Flowever, in spite of this 
spatial mismatch between the model single point simulation 
and available measurements in multiple locations around the 
model site, with few exceptions, the overall agreement be- 
tween model and observations is quite good. It will be shown 
later in this section that the majority of model versus ob- 
served parameters are statistically correlated with /7-values 
well below 0.05. 

A quantitative approach based on metrics, assessment in- 
dices and skill scores is provided to evaluate model perfor- 
mance. For this purpose, we use a Taylor diagram (Jolliff et 
al., 2009) to conduct the model evaluation based on all avail- 
able in situ and satellite observations. The normalized stan- 
dard deviation (er*) and the correlation coefficient ( R) from 
the model (m) to reference field (r) comparisons may be dis- 
played on a single Taylor diagram (Fig. 6). Note that er* is 
displayed in both x and y coordinates of the polar diagram 
of Fig. 6. The Taylor diagram is a polar coordinate diagram 
that assigns the angular position to the inverse cosine of R 
(cos -1 (/?)). A correlation of zero is thus 90° away from a 
correlation of 1. The radial (along-axis) distance from the 


origin is assigned to the normalized standard deviation. The 
reference field point (black circle in Fig. 6), which is com- 
prised of the statistics generated from a redundant reference 
to reference comparison, is indicated for the polar coordi- 
nate (1.0, 0.0). The model to reference comparison points 
may then be gauged by how close they fall to the reference 
point. This distance is proportional to the normalized unbi- 
ased root-mean-square difference, as defined by the equation 
RMSD*’ = (I.0 + cr* 2 -2cr*R) 0 ' 5 , where a* = a m !a r . 

Inspection of Fig. 6 reveals that the great majority of 
model versus reference comparisons fall within significant 
R values (0.7-0. 9) and a* close to 1. Except for alkalin- 
ity (TA), all other R -values are significant at the 95 % con- 
fidence level, with p-values well below 0.05. The TA refer- 
ence data have a relatively low correlation with the model, 
probably as a result of the low local TA variability, even 
though the bias between model and observations is small. 
The Taylor diagram does not provide an estimate for the bias, 
but Table 2 provides surface ocean means and standard devi- 
ations (STDs) for all the variables shown in the Taylor dia- 
gram. The means and STDs for the model variables are very 
close to the observed values. 
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Fig. 8. Time series of model versus satellite-derived Chi a, PIC (b), and PP (c). The black lines are the daily model parameters. The black 
dots in (a) and (b) are monthly SeaWiFS standard products, while the black dots in (c) are monthly satellite-derived CbPM PP values. The 
blue line time series in panel (b) is the model-derived euphoric depth to mixed layer thickness ratio (Zeu/MLD), and the two green triangles 
are PP values from in situ C-14 incubations. Note that the y-scales for Chi and PP are logarithmic. 

Table 1. Tabulated values show comparisons between observed and modeled mean DIC, TA, 1026 and surface ocean pC0 2 . Values are for 
the reference (Experiment A), no relaxation (Experiment B), and abiotic (Experiment C) runs. Only daily model values matching the nearest 
observed dates were included. Differences between model experiments and observed values (Ad) are also included. The time periods for 
averaging are 1991-2007 for DIC, 1997-2007 for TA, and 1993-2007 for pC02. 


Variable 

Observed 

A 

Ad 

B 

Ad 

C 

Ad 

DIC 

2111.01 

2106.87 

-4.14 

2095.93 

-15.08 

2122.74 

+11.73 

TA 

2309.26 

2307.82 

-1.44 

2303.82 

-5.44 

2303.89 

-5.37 

pC0 2 

347.65 

350.61 

+2.96 

332.30 

-15.35 

391.55 

+43.90 


Despite having significant R values, the model has some- 
what larger nitrate er* than the observed reference. Time se- 
ries of observed and model nutrients (Fig. 9 in Sect. 6) re- 
veal that the model underestimates all nutrients during 2005- 
2006, especially nitrate. Flowever, the low winter nutrients 
are consistent with warmer winter SSTs, as given by both 
model and observations (Fig. 10 in Sect. 6), and consequently 
shallower MLDs. Since these are the only two years during 
which the model underestimates nutrients, it is possible that 
the monthly mean physical fields and monthly climatological 
nutrients used to derive the forcing for the model were of in- 
sufficient resolution for capturing the smaller scale variabil- 
ity responsible for nutrient enhancement that occurred during 
those two years. Surface ocean /?C0 2 references are in good 


agreement with the model ( R ~0.7 and er*~1.0). Finally, 
the model PP provides good agreement with the PP estimates 
derived with satellite -based CbPM (Fig. 8 in Sect. 6) using 
SeaWiFS data at the precise location of the model simulation 
{r 2 =0.68). 

5 Sensitivity runs 

In order to evaluate the effects of model forcing and biologi- 
cal uptake of carbon on the model results, a couple of sensi- 
tivity runs (Experiments B and C) were performed and com- 
pared with the baseline run (reference, Experiment A), which 
was conducted following the description given in 2.2. One of 
the sensitivity runs (Experiment B) was conducted without 
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Fig. 9. Time series of model surface nutrients (NO3, PO 4 , SiC>2) with corresponding in situ observations from Suratlante (red dots) and 
CARINA (blue dots). 


the relaxation of DIC and TA to the empirical values derived 
from the data. In addition, the deep water trend in DIC was 
turned off and deep water values of DIC and TA were main- 
tained constant and equal to their initial values. The second 
sensitivity run (Experiment C) consisted of maintaining the 
DIC and TA forcing strategy used in the reference run but 
with the biological effects turned off, i.e., abiotic ocean case. 

Four surface ocean carbon-related variables were ana- 
lyzed: DIC, TA, pCC> 2 , and air-sea CCb flux. The clima- 
tologic (1981-2008) seasonal variability of these variables is 
shown in Fig. 7a. There are notable differences between the 
three different model experiments. Although Experiment B 
(without relaxation and DIC bottom trend) differs slightly 


from the reference, the largest differences are apparent in 
the abiotic run (Experiment C), since the biological sources 
and sinks that affect the DIC and TA seasonality were turned 
off in the model. These are the net community production, 
which is a sink for DIC (CO 2 uptake) and source for TA (ni- 
trate uptake consumes hydrogen ions), dissolution (source 
for both DIC and TA) and CaCC >3 production (sink for both 
DIC and TA). With these terms turned off, only thermody- 
namic and vertical mixing effects are at work in the model 
and that is reflected in the seasonal changes of DIC, TA, 
pC0 2 , and CO 2 flux. Note that without the biology (Experi- 
ment C), the ocean becomes a source of CO 2 for most of the 
year, while the reference run (Experiment A), which agrees 
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Fig. 10. Time series of model SST and surface ocean DIC and pCC>2 . Corresponding in situ measurements are shown with their diverse 
origin color-coded. The surface ocean pC Cb (bottom tier) also shows the GLOBALVIEW atmospheric pCC>2 and the ocean pCC>2, derived 
from CARINA TA, DIC, temperature and salinity. 


best with observations (see next section), shows an ocean 
that is a sink of CCb for most of the year and a small source 
during April-May. This result agrees with previous studies 
that show that biology is critical to these seasonal changes 
at this location and that vertical mixing also plays a major 
role (Takahashi et al., 1993, 2002, 2009; Ullman et al., 2009; 
Bennington et al., 2009; Metzl et al., 2010). The annual mean 
values of / 2 CO 2 and air-sea CO 2 flux derived from the ref- 
erence, no relaxation, and abiotic runs (Experiments A, B, 
and C, respectively) are 337.0, 325.6, 371.0 patm, and +7.3, 
+12.5, —1.7 mmolCm -2 d -1 , respectively. Further compar- 
ison between mean DIC, TA, and surface ocean pC0 2 from 
observed values and model experiments shows significant 


differences between Experiments A, B and C when com- 
pared to observations. In this comparison only daily values 
nearest to the time of observations were included. The time 
periods for averaging were 1991-2007 for DIC, 1997-2007 
for TA, and 1993-2007 for pCC> 2 . The results are summa- 
rized in Table 1, which shows significant increase in the bi- 
ases for the no relaxation (Experiment B) and abiotic (Exper- 
iment C) runs when compared to the reference run (Experi- 
ment A). For example. The biases for surface ocean pCOj 
are +3 pa tin. —15 patm, and +44 patm for Experiments A, B, 
and C, respectively. 

The difference in interannual variability of surface DIC 
and ocean pCOi from two of the model experiments, one 
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Table 2. Summary of surface-ocean observed and model (reference run) mean values and standard deviations (STD) for a selected number 
of variables. Only daily model values matching the nearest observed dates were included. 


Variable 

Observed 
Mean STD 

Model 

Mean STD 

N 

SST (°C) 

9.02 

1.77 

8.92 

1.71 

372 

TA (//mol kg -1 ) 

2309.26 

8.88 

2307.82 

2.61 

185 

DIC (//mol kg - 1 ) 

2111.01 

24.42 

2106.87 

27.53 

243 

N0 3 (//M) 

9.38 

3.27 

7.00 

4.68 

102 

Si0 2 (//M) 

3.36 

2.14 

3.82 

2.26 

103 

P0 4 (gM) 

0.66 

0.21 

0.62 

0.22 

122 

pC0 2 (//atm) 

347.65 

26.21 

350.61 

31.18 

199 

Sat PP (mgC m -2 d -1 ) 

207.85 

323.19 

220.63 

235.28 

114 



Fig. 11. Seasonal cycles of in situ CPR cell counts for diatoms, 
dinoflagellates, and coccolithophores (a), and model nitrogen-based 
concentrations of the same functional groups (b). 


including deep water DIC and TA relaxation and bottom DIC 
trend (Experiment A) and one without it (Experiment B), is 
shown in Fig. 7b. The available data are also shown for com- 
parison. There is a striking difference between these two ex- 
periments. Experiment A shows a very good match with the 
available data, while the results from Experiment B have a 
very poor agreement with the data, including a significant re- 
duction in the seasonal amplitude of DIC and pCOj ■ This 
result illustrates the degree to which the model results are di- 



Fig. 12. Seasonal cycles of Chi a and calcite (a), PAR and SST (b), 
Zeu:MLD ratio (c), and nutrients (d). 
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Fig. 13. Seasonal concentration of phytoplankton functional groups 
biomass with and without the presence of coccolithophores (a and 
b, respectively), and net community production with and without 
coccolithophores (c). 

rectly a function of the setting of DIC, TA and nutrients at 
>200 m, how much the deep reservoir impacts surface be- 
havior, and more clearly that it is very important to be in- 
cluded in the model to appropriately represent the interan- 
nual variability driven by non-local processes. 

6 Discussion 

6.1 Seasonal and interannual variability 

Figure 8 shows a 1998-2008 time series of model and 
satellite-derived Chi a, calcite concentration (PIC), and PP. 
The ratio of euphotic depth to mixed layer depth (Zeu/MLD) 
derived from the model is shown by the blue line, super- 
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posed to the PIC time series. This ratio provides a mea- 
surement of how much of the daily light exposure the phy- 
toplankton experiences during the growth season. Appar- 
ently, every year the most productive periods (when Zeu is 
shallowest) occur during highly stratified conditions (with 
shallow MLDs). That is an indication that light availabil- 
ity for photosynthesis is a major limiting factor, causing sea- 
sonal variability of phytoplankton growth and probably also 
affecting the species succession in the region (Nanninga and 
Tyrrell, 1996; Tyrrell and Taylor, 1996). 

Overall, the model agrees well with the corresponding 
satellite-derived products. In fact, the model reproduces the 
satellite PIC interannual variability rather well, albeit with 
the model showing a tendency to produce peak PIC concen- 
tration higher and later than the satellite observations. No- 
tably, the interannual changes in peak PIC concentrations are 
correlated with light intensity, as they follow the variability 
of the Zeu/MLD ratio. The variability in this ratio is pri- 
marily driven by changes in MLD, since the standard devia- 
tion in MLD during the summer months (June-July-August, 
1981-2008) is about 10 times larger than the standard devi- 
ation in Zeu (24.8 ± 10.8 m and 72.6 ± 1.4 m, respectively). 
The growth of calcite-forming coccolithophores in the model 
is, by design, three times more demanding on light exposure 
than the other two functional groups, since it has the low- 
est assigned P vs. I initial slope of 0.033 (Wm -2 ) -1 day -1 
(see Table B2). The coherent model prediction of growth of 
coccolithophores and the variations in PIC concentration are 
consistent with the interannual changes in PIC observed by 
the satellite, which confirms the relevance of the Zeu/MLD 
ratio. 

Time series of model and in situ measurements of nutri- 
ents (NO3, PO4, and Si0 2 ) for the time period, for which 
there are data (1990-2008), are shown in Fig. 9; SST and 
surface ocean DIC and pCOi for the same time period are 
shown in Figure 10, also with corresponding in situ mea- 
surements. The GLOBALVIEW atmospheric pCOj, which 
was used to calculate A/jCCH and sea-air CO2 flux, is shown 
in the bottom panel of Fig. 10 along with the surface ocean 
pC0 2 . The seasonal and interannual changes in nutrients and 
carbon are substantial and driven primarily by changes in bi- 
ological production, vertical mixing, SST, and sea-air CO2 
flux. 

The observed and model seasonal partitioning of phy- 
toplankton functional groups (diatoms, dinoflagellates, and 
coccolithophores) are shown in Fig. 11. The observed sea- 
sonality originates from the aggregated, monthly-averaged 
indices of phytoplankton groups from the CPR data (cell 
counts) in the standard area B6, south of Iceland (Fig. 11a). 
The size of standard area B6 is relatively large (59° N to 
64° N; 19° W to 31° W), so some spatial smoothing is to be 
expected. The model shows (Fig. 11a) that diatoms are the 
most abundant functional group, with a large bloom in early 
spring and a secondary bloom in fall. The diatom bloom is 
followed by blooms of dinoflagellates and coccolithophores. 
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Fig. 14. Seasonal (1998 — 2008) changes of alkalinity (a), DIC (b), surface ocean pCC >2 (c), and sea-air CO 2 flux (d) with (red) and 
without (blue) coccolithophores. The black line is the difference. 


The double diatom peak is also present in the model nitrogen- 
based concentrations of the functional groups (Fig. 1 lb). The 
abundance of coccolithophores in the CPR data is only a 
fraction of what one may expect, due to the small size of 
the organisms (5 pm) compared to the size of the mesh used 
(250 pm). The CPR data are however a relevant index of 
the relative annual changes and thus reveal the blooming of 
coccolithophores. The nitrogen-based concentration of coc- 
colithophores in the model is thus much more substantial and 
in line with the findings of previous field work studies in the 
northeast North Atlantic, reporting large concentrations of 
coccolithophores during the summer bloom. For example, 
Fernandez, et al. (1993) reported a total estimated coccol- 
ithophore C biomass as large as 50 % of the total phytoplank- 
ton C biomass. 

6.2 Biogeochemical response to physical forcing 

The seasonal changes in biogeochemical properties as a re- 
sult of physical forcing are given in Fig. 12, which shows 
the model seasonal cycles of chlorophyll, calcite (PIC), SST, 
PAR, the ratio of euphotic depth to MLD (Zeu/MLD), and 
nutrients (NO3, PO4, and SiCb). The rise of Chi a starts in 
April, when the Zeu/MLD ratio is ~0.5 and a drawdown of 
nutrients ensues. That is as one may expect, according to the 
Sverdrup (1953) theory on vernal blooming of phytoplank- 
ton, as the critical depths in oceanic environments are several 
times that of the Zeu. 

The PIC concentration starts to rise in May when the 
Zeu/MLD has reached ~2. The broad peak in the total 
Chi a from all three functional groups extends from May 
to September, associated with shallower MLDs (Zeu/MLD 
ratio > 2.0). Depletion of available light and nutrients for 
growth and an increase of grazing pressure reduce biomass 


significantly during the summer. The growth season cycle is 
completed during autumn, as wind stress increases the verti- 
cal mixing, promoting nutrient renewal. The coccolithophore 
bloom, indicated by the changes in PIC concentration due to 
calcite production, peaks in July-August, when light condi- 
tions (Zeu/MLD ratio > 4.0) provide favorable conditions 
for the bloom-forming coccolithophore Emiliania huxleyi, a 
species well-known to the area (Balch et al., 1992; Fernandez 
et al., 1993; Holligan et al., 1993a; Holligan et al., 1993b). 
As the MLD deepens and light (PAR) levels are significantly 
reduced, Chi a and PIC concentrations drop gradually after 
September. The model confirms established ideas of the in- 
terplay of light availability, vigorous winter mixing/summer 
restratification, and nutrient availability, which are central 
to the classical North Atlantic spring bloom (Nanninga and 
Tyrrell, 1996; Tyrrell and Taylor, 1996). 

6.3 Impact of phytoplankton blooms on carbon uptake 

To evaluate the impact of phytoplankton groups on the up- 
take of carbon, we conducted a model experiment with the 
coccolithophore components turned off (no carbonate pump) 
and compared the results with the baseline experiment (bio- 
logical and carbonate pumps), including all three functional 
groups (diatoms, dinoflagellates, and coccolithophores). The 
results are shown in Figures 13 and 14 as seasonal 8-day cli- 
matologic averages for 1998-2008. In Figure 13 we see the 
results of seasonal changes in the upper 120 m vertically- 
integrated phytoplankton biomass concentrations with and 
without the presence of coccolithophores (Fig. 13a and b, 
respectively), and the corresponding changes in net com- 
munity production (Fig. 13c). In the presence of coc- 
colithophores, the yearly-averaged integrated biomass for 
diatoms, dinoflagellates, and coccolithophores was 534, 
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318, and 184 mg Cm -2 , respectively. The correspond- 
ing model-derived Np with all three functional groups was 
134 mg C m -2 d - 1 . With the exclusion of coccolithophores, 
the population of diatoms increased to an integrated value 
of 613 mg Cm -2 and dinoflagellates to 341 mg Cm -2 , re- 
spectively. However, there was a net decrease in Np to 
104 mg C m -2 d - 1 (or a 22.4 % decrease). At first glance, it 
would be expected that the diatoms and dinoflagellates would 
benefit from the absence of coccolithophores and have ad- 
ditional nutrients to draw from and increase their combined 
production. However, the explanation lies in the carbon up- 
take efficiency of coccolithophores when compared to the 
other two groups. In equation B35 of Appendix B, it can 
be seen that, for the baseline run, dinoflagellates and diatoms 
have a Redfield C:N (6.625), while coccolithophores have a 
higher C:N (9.4), which is an average, based on the range 
of 5.81 to 13.05 reported by Fernandez et al. (1993). The 
diatoms and dinoflagellates do pick up some of the slack 
but not entirely, because the coccolithophores uptake carbon 
42% more efficiently ( 1 00[ 1 -9.4/6.625)). In addition, the 
peak demand of nutrients does not overlap among the three 
functional groups, as the maximum productive time periods 
of each group are not in phase. So the absence of coccol- 
ithophores is less influential on the availability of nutrients 
for the other two groups. This seems to be a consistent ex- 
planation of why the removal of coccolithophores from the 
model results in the reported differences in PP and pC0 2 
drawdown. 

A couple of sensitivity tests were conducted to assess 
the effect of changing the C:N ratio for coccolithophores 
in the model. We used two values of C:N ratio for coccol- 
ithophores in the model, Redfield (R1 = 6.625), and twice 
the Redfield ratio (R2 — 13.25). The annual mean (1998 — 
2008) results from R1 for the following key model pa- 
rameters, pC0 2 , DIC, TA, air-sea CCD flux, net commu- 
nity production, CaCC >3 concentration, and PIC produc- 
tion are 358.8patm, 21 14.6 pmol kg -1 , 2308.3 pmolkg -1 , 
4.4 mmol Cm -2 d -1 , 123.6 mg C m -2 d -1 , 3.3 mgm -3 , and 
4.9mgCm -2 d -1 . The equivalent values in the same units 
for experiment R2 are 353.7 (—5.1), 2111.3 (—3.3), 2308.0 
(-0.3), 5.73 (+131%), 146.1 (+118%), 6.7 (+200%), and 
9.8 (+200%). So the model is very sensitive to changes in 
the C:N ratio for coccolithophores. 

The 22% reduction in Np between the baseline run 
(C:N = 9.4 for coccolithophores) and the baseline run without 
coccolithophores has consequences in the biological uptake 
of carbon. Figure 13 shows that, without the net commu- 
nity production and calcification by coccolithophores, at the 
peak of the biological drawdown, the alkalinity (Fig. 14a) in- 
creases by up to 4 pmol kg - 1 with a corresponding increase 
in DIC (Fig. 14b) of up to 16 pmolkg -1 . As a result, the net 
effect of the absence of the coccolithophores bloom is an in- 
crease in pC0 2 of more than 20 patm (Fig. 14c). The largest 
differences occur from July to October, which translate into 
ApCQ 2 differences ranging from 24% to 35% during this 
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period. The seasonal changes in the sea-air COt flux are 
shown in Fig. 14d. At the peak of the bloom, the exclusion of 
coccolithophores causes a reduction of CO 2 uptake of more 
than 6 mmol m -2 d -1 . On average, the impact of coccol- 
ithophores on air-sea CO 2 flux is 27 %. This result highlights 
the importance of including all major functional groups in 
the modeling of carbon variability in the subpolar North At- 
lantic. Some sensitivity tests were conducted using different 
gas transfer coefficients (k). The baseline run uses a value of 
0.31 (Wanninkhof, 1992), as indicated in Eq. (B39) in Ap- 
pendix B. A value of k — 0.27 (Sweeney et al., 2009) was 
used to evaluate the model sensitivity of air-sea CO 2 flux to 
different gas transfer coefficients. A reduced CO 2 ocean up- 
take of ~7 % is obtained with the experiment using k = 0.27, 
when compared with the experiment using k=0.31, with both 
experiments including coccolithophores. Model experiments 
using k — 0.27 with and without coccolithophores indicate 
that the exclusion of coccolithophores causes a reduced CO 2 
ocean uptake of ~28%, a small difference (1 %) compared 
to the use of k = 0.31. 

Suykens et al. (2010) reported the carbon uptake impact 
of coccolithophore blooms in the northern Bay of Biscay 
based on DIC and alkalinity observations. They concluded 
that the decrease of DIC (and increase of pC0 2 ) due to net 
community calcification was overwhelmingly lower than the 
decrease of DIC (and decrease of pC0 2 ) due to net com- 
munity production, a clear indication of the importance of 
coccolithophore blooms in the uptake of atmospheric CO 2 in 
the North Atlantic, where blooms of coccolithophores are the 
most intense and recurrent. 

A spatial assessment on the reduction in uptake of atmo- 
spheric CO 2 under the exclusion of coccolithophores was 
made for the regional domain of the Irminger and Icelandic 
Basins. The region considered is bound by 44° W - 8° W 
and 56° N - 65° N, and by depths greater than 200 m to ex- 
clude coastal regions of Iceland and Greenland. The monthly 
areal (in m 2 ) extent of the bloom was determined using Sea- 
WiFS monthly composites within the above boundaries and 
considering only pixels with PIC concentrations larger than 
5 mgm -3 , a value arbitrarily picked as the threshold of sig- 
nificant influence of coccolithophore blooms in the draw- 
down of atmospheric CO 2 (monthly concentrations range 
from 0 to 27 mgm -3 ). For comparison, a threshold of 
2 mg m -3 in PIC concentration increases the flux estimate by 
about 12% (or about 1 80 tonnes Cyr -1 ). The monthly CO 2 
flux difference was derived from the product of satellite area 
and the flux difference shown in Fig. 14d. The annual mean 
flux difference was estimated at nearly 1500 tonnes C yr - 1 . 

The seasonal drawdown of surface ocean pC0 2 is a result 
of two competing effects, i.e., temperature warming and bi- 
ological uptake effects. Takahashi et al.(2002) developed a 
method to separate these two effects for the global oceans. 
The effect of biology ( B e ) on the surface ocean / 2 CO 2 in a 
given area is represented by the seasonal amplitude of pC0 2 , 
corrected to the mean annual temperature in that area. The 
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effect of temperature changes {T e ) on the seasonal pCOj 
variations is represented by the seasonal amplitude of the 
mean annual /;C0 2 , corrected to the range of observed tem- 
peratures. 

B e = (ApC02)bio = (pC0 2 at 7mean)max (P^O 2 at Tmeanlmin 

T e = (ApC0 2 ) temp = (pC0 2 atr obs ) max - (pC0 2 atr o b s ) min (1) 

where the subscripts “min” and “max” indicate the seasonal 
minimum and maximum values, 7 0 |, s is the observed temper- 
ature for each p CO2 value and r mea n is the mean tempera- 
ture over the entire seasonal cycle. The relative importance 
of the biology and temperature effects can then be expressed 
by the ratio B e /T e . Using climatologic (1981-2008) model 
seasonal surface ocean 72CO2 and SST, we calculate the ratio 
B e /T e as 1.92, which means that the biology effect is nearly 
twice the temperature effect in the shaping of the seasonal 
/;C0 2 variability at the model location. 

Although it is not straightforward to infer to what extent 
the model results from a single location in the subpolar North 
Atlantic are representative of the North Atlantic in general, 
previous studies have shown coherence of results in much 
broader regions of the North Atlantic. For example, McKin- 
ley et al. (2011) have shown that the North Atlantic can be 
separated in three large biogeographic regions (“biomes”) in 
terms of surface ocean p C0 2 variability. These biomes were 
assigned on the basis of annual maximum mixed layer depth, 
annual mean SeaWiFS Chi a, and SST. Namely, the North 
Atlantic was divided into a northern seasonally stratified gyre 
(SP-SS) biome, a southern permanently stratified subtropi- 
cal gyre (ST-PS) biome, and a seasonally stratified subtrop- 
ical (ST-SS) biome. For a 29yr-long observational period 
(1981-2009), the surface ocean /2CO2 trends converge with 
the atmospheric growth for all 3 biomes, which occupy 87 % 
of the total area of the North Atlantic. When a shorter time 
scale (1993-2005) is considered, the three biomes have dif- 
ferent trends, but the SP-SS biome itself has a significant ge- 
ographic extension, including a large region south of Iceland, 
and the Irminger, Icelandic, and Norwegian Seas. Our 1-D 
model location lies within the SP-SS biome and therefore is 
representative of surface ocean pC0 2 decadal variability for 
that region. Flowever, for shorter time scales, more signifi- 
cant regional differences are likely to occur. 

7 Summary and conclusions 

A 1-D ecosystem model was developed to assess the rela- 
tive contribution of phytoplankton functional types to atmo- 
spheric C0 2 uptake. Using a Taylor diagram, skill assess- 
ment of model versus field measurements reveals high scores 
for the majority of biogeochemical parameters, for which in 
situ data are available. 

The seasonal patterns of phytoplankton concentrations are 
a response to the interplay between light availability, vig- 
orous winter mixing/ summer restratification, and nutrient 


S. R. Signorini et al.: The role of phytoplankton dynamics 

availability, not unlike the classical North Atlantic spring 
bloom (Ducklow and Flarris, 1993; Weeks et al., 1993). 
Functional groups compete seasonally for ideal growth con- 
ditions (Litchman et al., 2007; Sambrotto et al., 1993; Sier- 
acki et al., 1993). Model results indicate that the spring- 
summer bloom consists predominantly of diatoms, with still 
significant but less intense blooms of dinoflagellates and coc- 
colithophores. The model shows that the diatom biomass 
peaks in May, with a secondary and less intense bloom in 
September. The dinoflagellates and coccolithophores peak 
in July through August, during which drawdown of surface- 
ocean C0 2 reaches its maximum value. The effect of biolog- 
ical changes in the surface ocean pC0 2 exceeds the temper- 
ature effect by a factor of almost 2, a clear indication of the 
importance of phytoplankton photosynthesis on the uptake of 
atmospheric C0 2 in the region, a result that is in agreement 
with previous studies (Takahashi et al., 1993, 2002, 2009; 
Ullman et al., 2009; Bennington et al., 2009; Metzl et al., 
2010 ). 

Model experiments were conducted to investigate the sea- 
sonal changes in phytoplankton concentration with and with- 
out the presence of coccolithophores, and their impact on 
carbon uptake. Without the influence of coccolithophore 
blooms, the alkalinity increases by almost 4pmolkg -1 and 
DIC is elevated by up to lhpmolkg^ 1 . The net effect of 
coccolithophores blooms is a drawdown in pC 02 of up to 
about 20 pa tm during summer, with a corresponding increase 
of atmospheric C0 2 uptake of about 6 mmolm -2 d -1 , an in- 
dication of the importance of including all major phytoplank- 
ton functional groups when modeling the biological carbon 
pump variability in the subpolar North Atlantic. 

Appendix A 

Ice-ocean model description (GPOMZ) 

The ocean model is hydrostatic and Boussinesq and uses a 
generalized vertical coordinate system, as described by Mel- 
lor et al. (2002) with a modified scalar advection scheme to 
avoid overshooting at sharp fronts (Mauritzen and Hakkinen, 
1997). The equation of state is formulated in terms of in situ 
density (Mellor, 1991), expressed as a function of potential 
temperature, salinity, and pressure. The model’s prognos- 
tic variables are the horizontal velocity components, poten- 
tial temperature and salinity, and twice the kinetic energy, q 2 
and q 2 l, where l is the turbulence macroscale (Mellor and 
Yamada, 1974, 1982). These turbulence quantities, together 
with the vertical velocity shear and buoyancy, determine the 
vertical mixing coefficients for momentum and scalar vari- 
ables. The dynamic-thermodynamic ice model is coupled to 
the ocean model via interfacial stresses and via salinity and 
heat fluxes through the ice -water interface. The thermody- 
namic component of the ice model follows Semtner’s (1976) 
formulation but with modifications to account for leads. The 
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Fig. Al. Maps of GPOMZ (a) and MODIS Aqua (b) annual mean SST for 2009. Figure Ala shows the geographic domain of GPOMZ, with 
two open boundaries: one across the subtropical South Atlantic at 16.7° S and one across the Bering Strait in the Pacific sector of the Arctic 
Ocean. The satellite SST in Fig. Alb is included for comparison, showing that the model provides a realistic representation of all surface 
fronts. Regions with SST less than zero are partially to fully ice covered. 



Fig. Al. Time series of monthly GPOMZ (black) and MODIS 
(blue) SST at the location of ECO1D-2.0 simulation (60° N, 
30° W). GPOMZ temperature profiles were bias-corrected based on 
monthly WOA05 climatologic profiles. GPOMZ SST is 0.78 °C 
warmer on average (1981-2008) when compared to WOA05 SST 
climatology. The adjusted GPOMZ SST has a 0.31 °C warm bias 
when compared to MODIS SST (July 2002 - December 2008), as 
shown in this figure. 


dynamics of the ice model is described as a continuum with 
a generalized viscous rheology (Hakkinen, 1987; Hakkinen 
and Me 11 or, 1992). 

The heat and momentum exchanges between the atmo- 
sphere and ice/ocean are computed from bulk formulae using 
an exchange coefficient of 1.4 x 10~°, whereas the momen- 
tum, salt, and heat exchanges between the ice-ocean inter- 
face and the ocean interior are computed from the “law of 
the wall”, assuming a turbulent boundary layer over a rough 
surface. The scheme uses the Yaglom and Kader (1974) 
formulation for heat and mass transfer over a rough wall. 



Fig. A2. Seasonal climatology of mixed layer depth (MLD) from 
de Boyer Montegut et al. (2004) (blue) and from GPOMZ at the 
ECO1D-2.0 site. There is good agreement from May to September, 
but GPOMZ MLD is much deeper in winter and late fall. 

which gives different salt and heat resistance coefficients and 
a Smagorinsky scheme for the lateral viscosity and diffusiv- 
ity . The vertical diffusivities for momentum and scalar quan- 
tities are computed from the turbulence closure and added to 
a background diffusivity of 1.0 x 10 -5 m 2 s _1 . The Jerlov 
formulation for oceanic optical water type IA (Jerlov, 1976) 
is used wherein 3 1 % of the incoming short-wave radiation 
is absorbed below the sea surface, decreasing exponentially 
with depth with an extinction coefficient of 0.042 m -1 . The 
model equations are solved using a split mode scheme and 
a z-level coordinate, curvilinear grid. The coupled ice-ocean 
model has a variable resolution of 28-56 km in an orthogo- 
nal curvilinear grid designed for the area. This resolution is 
not eddy-resolving, but is adequate for simulating large-scale 
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circulation. The model domain boundary in the Atlantic sec- 
tor lies at 16.7° S and at the Pacific side in the Bering Strait. 
The model topography is a smoothed version of the ET0P05 
digital database (ET0P05, 1988). The model uses a total of 
26 vertical levels as follows: 3, 9, 18.5, 37.5, 75, 125, 175, 
250, 350, 450, 575, 725, 900, 1100, 1300, 1500, 1700, 1900, 
2150, 2450, 2800, 3250, 3750, 4250, 4750, and 5250 m. Fig- 
ure Al shows maps of GPOMZ and MODIS Aqua annual 
mean SST for 2009. Figure Ala shows the geographic do- 
main covered by GPOMZ. The satellite SST in Fig. Alb is 
included for comparison, showing that the model provides 
a realistic representation of all surface fronts. Regions with 
SST less than zero are partially to fully ice covered. Fig- 
ure A2 shows a time series comparison between GPOMZ and 
MODIS SST at the location of ECO1D-2.0 simulation. Fig- 
ure A3 shows the monthly seasonal MLD climatology from 
de Boyer Montegut et al. (2004) and from GPOMZ at the 
ECO ID location. There is good agreement from May to 
September, but the GPOMZ MLD is much deeper for the 
other months. The deBoyer-Montegut MLD is provided on 
a 2° x 2° resolution grid with three different MLD criteria. 
The criterion for both MLDs in Fig. A3 is the fixed density 
criterion, i.e., the depth where density increase compared to 
density at 10m equals 0.03 kgm -3 . 

The ocean model was initialized with the winter average 
hydrographic climatology derived from the Polar Science 
Center Hydrographic Climatology (PHC 3.0) (Steele et al., 
2001). A transport of 0.8 Sv was specified at Bering Strait 
and 0.8 Sv out at the southern boundary. At the northern 
and southern boundary, the salinities and temperatures are 
relaxed to monthly climatological values. Climatological 
monthly river discharge for 474 rivers within the model do- 
main is specified with the global river database, described by 
Barron and Smedstad (2002). An additional river discharge 
of 700 km 3 yr _1 is uniformly distributed along the Eurasian 
margin of the Arctic. River runoff is treated as a virtual flux 
of salt. The model is forced with daily NCEP/NCAR Reanal- 
ysis sea level pressure, air temperature 2 m above sea level, 
surface momentum fluxes and relative humidity for 1 948— 
2011. Climatological monthly mean cloudiness is based on 
the ISCCP D2 data set from 1984-2004, except north of 
~65° N where a uniform monthly mean cloudiness is speci- 
fied. Climatological monthly mean precipitation minus evap- 
oration (P— E) is derived from the NCEP operational analysis 
(Rasmusson and Mo, 1996). An additional 1350km 3 yr _1 of 
precipitation is added, and distributed with a peak value at 
~8° N to simulate the Intertropical Convergence Zone. 


(N0 3 , P0 4 , NH 4 , Si0 2 , and Fe), POC, DIC, DOC, alkalin- 
ity (TA), calcite production, chlorophyll, complete carbon- 
ate chemistry and air-sea C0 2 flux. The governing equa- 
tions for the biogeochemical model are provided hereafter, 
where the subscripted index i = 1,2,3 represents diatoms, 
dinoflagellates, and coccolithophores, respectively. Iron (Fe) 
limitation, although included in the model, was not consid- 
ered in this study However, there is evidence (Nielsdottir et 
al., 2009) of iron limitation of the post-bloom (July to early 
September) phytoplankton communities in the Iceland Basin, 
east of the model location, where high nutrient-low chloro- 
phyll (HNLC) conditions may occur. Iron limitation studies 
in the Iceland Basin will be a topic of future studies using the 
same model. 

The basic form of the model equations for a generic tracer 
of concentration C, differentiated with respect to time t and 
vertical coordinate z and balanced by the sources and sinks 
of the particular biogeochemical property, is: 


dC dC 9 T 3 C' 

h w k v 

3 t dz 3 z 3 z 


sources - sinks 


where w is the vertical velocity and k v is the vertical eddy 
diffusivity. We simplify the notation of the left hand side of 
all equations by substituting it with the single term 


B1 Phytoplankton (diatoms, dinoflagellates, 
and coccolithophores) 


d Pi ^ 

^ — Mi.phy-* i Mi, zoo Z (^phy,N0 3 ”1“ ^phy,det) 

3 Pi 

(Pi ~ Pi.o ) — WPi 

oz 


(Bl) 


B2 Zooplankton 

— j — y y ' (/n.zooZ) — ^zoo.M^tz — z 0 ) M ZOOf( \ e i (Z — Z Q )~ (B2) 


B3 Nutrients 


dN0 3 1 

— -- -(N0 3 *-N0 3 ) 
d/ r 

— 2 ^ ^ (/A, phy Pj ) T~ A h 


nh 4 

0.07 + NH 4 


(B3) 


Appendix B 

Biogeochemical model description (ECO1D-2.0) 

The model features multiple functional groups (diatoms, di- 
noflagellates, and coccolithophores), zooplankton, nutrients 


dNH 4 x - 

— CIpM phy.det / . (Pi Pi,o) 4“ i^phy.NO. (B4) 

i 

^ ' (F*i Pi, o) ^ (Mi.phy-^z) remdet/y + Cl z M zoo,detZ 

i i 


H“^zoo,N 03 (Z Z 0 ) A n 


( NH 4 \ 
V0.07+NH4 / 


-b k m DON 
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Table Bl. Ecosystem/carbon model state variables definition, forcing strategy, and initial values. 


Symbol 

Units 

Parameter 

Forcing/Initial value 

no 3 

mmol N m -3 

Dissolved inorganic nitrate 

**Deep water NR 

P 0 4 

mmol Pm -3 

Dissolved inorganic phosphorus 

**Deep water NR 

Si 0 2 

mmol Si m -3 

Dissolved inorganic silica 

**Deep water NR 

nh 4 

mmol N m -3 

Ammonium 

0.05 

0 2 

/iiuol kg - 1 

Dissolved oxygen 

268.0 

CaC 0 3 

mmol C m -3 

Calcium carbonate (calcite) 

0.01 

DIC 

/Ltmol kg - 1 

Dissolved inorganic carbon 

**Deep water NR+trend 

TA 

/imol kg - 1 

Alkalinity 

**Deep water NR 

P, 

mmol N m -3 

Phytoplankton (( = 1 , 2 , 3 )* 

0.1 

Pi,o 

mmol N m -3 

P background concentration 

0.001 

Z 

mmol N m -3 

Zooplankton 

0.1 

Zo 

mmol N m -3 

Z background concentration 

0.001 

DON 

mmol N m -3 

Dissolved organic nitrogen 

1. 0/0.0 

DOP 

mmol Pm -3 

Dissolved organic phosphorus 

0 . 06 / 0.0 

DOC 

/iiuol C kg -1 

Labile dissolved organic carbon 

15 . 0 / 0.0 

det/v 

mmol N m -3 

Detrital particulate nitrogen 

0.001 

detp 

mmol Pm -3 

Detrital particulate phosphate 

0.001 

det 5 

mmol Si m -3 

Detrital particulate silicate 

0.001 

detpcaC 03 

mmol CaC 0 3 m -3 

Detrital particulate calcite 

0.001 


*Diatoms (i = 1), dinoflagellates (i — 2), and coccolithophores (i — 3) **Deep layer Newtonian relaxation to values in (B40), (B41), and (B42). 


dP 0 4 1 

— (P 0 4 * — PO4) = ^rp DOP + rem detp 

at t 

+[(fl/ 7 Afphy 5 det H - ^phy,N 03 ) I>-^> 

i 

(t^i,phyP0 “P &zMzoo,det Z £zoo,NC>3 (Z Z 0 )\ (B5) 

i r ^ 

- I (Si 0 2 * - Si 0 2 ) (B6) 

— Afphy 5 det "P ^phyjNC^X^l P\,o) 

^ ' (M2,phyp2) | _- l - ^z^zoo,det^ H“ ^zoo,N 03 (Z Z 0 ~) 1 ( rr j 
/ L J V ^ / Red 

-\-k rp DOS + rem dets 

B4 Dissolved organic matter 

= ftp Mphy,det ^ ' (P; Pi,o) ~ a z M Z oo,detZ k rn DON (B7) 

i 

^ = —[ftp Mphy,det } ' (Pi ~ Pi,o) + a z ^zoo,detZ^] — &rp DOP (B8) 

i 

“T = — [0-15 ^~^ (f J 'i.ph\Pi') ~^ a p ^phy.det ) ' (Pi ~ Pi,o)~\~Qz ^zoo.defZ 2 ] (B9) 

P i i 

— k TC DOC 

V N / Red 


B5 Detritus 
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+ (1 - a z - al)M zooAet Z 2 - Wdet^f 11 - rem det N 

dd*r = (1 - ap -flp[M phy , de tE^/ - Pi.o)-^] 

i 

+[(1 - y) Y] (Vi,zooZ) + (1 -a z - a')M z 00 , det Z 2 ] — 

" r 1 

3detp 
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3z 

^ = [(1 - a p - apMphy.detCA - P Uo ) 
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ddetcaeo 3 = |p _ 0p _ a pM phy , det(P 3 - P 3i0 ) 

"t“(l y )(/t- 3 ,zooZ) +(1 <U)/V/ ZO o,detZ~J 

3 det Ca co 3 

-«'det 

3 z 

B6 Calcite 

dCaC0 3 * 9CaC0 3 

~r = / J CaC0 3 - A. CaC03 CaC0 3 - W’det — 
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(B 14 ) 


www.geosci-model-dev.net/5/683/2012/ 


Geosci. Model Dev., 5, 683-707, 2012 



702 


S. R. Signorini et al.: The role of phytoplankton dynamics 


Table B2. Summary of ECO1D-2.0 parameters used in the 1981-2008 simulation for the subpolar North Atlantic. Equivalent parameters 
used in ECO ID- 1.0 for the BATS site are included in parenthesis where appropriate (ECO ID- 1.0 combines all phytoplankton into one single 
group, so phytoplankton- relevant parameters are provided under the diatoms group in this table). 


Parameter 

Symbol 

Units 

Value 

Diatoms (Pi) 




Maximum growth rate 

Pphy,max 

day -1 

1.2 (3.0) 

Initial P-I slope 

a 1 

(Wm -2 ) -1 day -1 

0.09 (0.25) 

Nitrate half saturation constant 

^N03 

mmol N m -3 

2.0 (0.025) 

Phosphate half saturation constant 

hoA 

mmol P m -3 

0.0015 (0.0025) 

Silicate half saturation constant 

%02 

mmol Si m -3 

0.3 

Ammonium half saturation (all functional types) 

*NH4 

mmol N m -3 

0.5 (0.025) 

Respiration rate 

R phy,N03 

day -1 

0.05 (0.05) 

Mortality rate 

^phy.det 

day -1 

0.05 (0.05) 

Maximum sinking speed 

Wp 

m day - 1 

1.0 

Dinoflagellates (P 2 ) 

Maximum growth rate 

Mphy,max 

day 1 

0.65 

Initial P-I slope 

«1 

(Wm -2 ) -1 day -1 

0.10 

Nitrate half saturation constant 

^N03 

mmol N m -3 

0.67 

Phosphate half saturation constant 

*P04 

mmol P m -3 

0.0015 

Respiration rate 

R phy,N03 

day -1 

0.05 

Mortality rate 

'^pliy.dct 

day -1 

0.05 

Coccolithophores (P 3 ) 

Maximum growth rate 

Mphy,max 

day 1 

1.15 

Initial P-I slope 

a 1 

(Wm -2 ) -1 day -1 

0.033 

Nitrate half saturation constant 

*-’N03 

mmol N m -3 

1.0 

Phosphate half saturation constant 

*P04 

mmol P m -3 

0.0015 

Respiration rate 

R phy,N03 

day -1 

0.05 

Mortality rate 

Mphy.det 

day -1 

0.05 

Maximum sinking speed 

Wp 

m day - 1 

5.0 

Zooplankton (Z) 

Maximum grazing rate 

M zoo, max 

day 1 

1.2 ( 1 . 0 ) 

Ingestion half saturation constant 

k phy 

(mmol N m -3 ) -2 

0.25 (0.25) 

Assimilation efficiency 

Y 


0.75 (0.75) 

Excretion rate 

^zoo,N03 

day -1 

0.1 (0.05) 

Mortality rate 

^zoo,det 

day -1 (mmol N m -3 ) -1 

0.1 (0.05) 


Dissolved organic matter (DOM), detritus (det), and remineralization 


Remineralization rate of detritus 

rem 

day -1 

0.01 ( 0 . 001 ) 

Remineralization rate of DON 

km 

day -1 

0.00001 ( 0 . 00001 ) 

Remineralization rate of DOP 

krp 

day -1 

0.00001 ( 0 . 00001 ) 

Remineralization rate of DOSi 

krs 

day -1 

0.00001 

Remineralization rate of DOC 

krc 

day -1 

0.0005 (0.0005) 

Fraction of dead phytoplankton converted to NH 4 

dp 


0.8 ( 0 . 8 ) 

Fraction of dead zooplankton converted to NH 4 

a z 


0.8 ( 0 . 8 ) 

Fraction of dead phytoplankton converted to DOM 

a 'p 


0.1 ( 0 . 1 ) 

Fraction of dead zooplankton converted to DOM 

< 


0.1 ( 0 . 1 ) 

Maximum sinking speed 

“Met 

m day - 1 

1.0 (2.5) 

Ammonium (NH 4 ) nitrification 

Maximum rate of ammonium nitrification 

^n,max 

// mol m 3 day 1 

0.02 ( 0 . 02 ) 

Minimum light inhibition dosage for nitrification 

^min 

Wm “ 2 

0.0095 (0.0095) 

Half saturation dosage for nitrification photoinhibition 

k d 

Wm ~ 2 

0.036 (0.8) 

Calcite (CaC 03 ) 

Maximum calcification to organic carbon production 

s CaC03 


0.4 

Calcite dissolution rate 

Y CaC03 

day 1 

0.03 

Chlorophyll a (Chi a) 

Minimum chlorophyll to nitrogen ratio 

Chl:N min 

mg Chi (mmol N) 1 

1.5 (0.8) 

Maximum chlorophyll to nitrogen ratio 

Chl:N max 

mg Chi (mmol N ) -1 

3.5 (1.6) 

Critical irradiance for photoadaptation 

/* 

Wm -2 

25.0 
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B7 Dissolved inorganic carbon 


dDIC 1 / \ 

— DIC — DIC ) = <5(z) 

d? r V / 

+ “ [^CaC 0 3 C a C 03 — -PcaC 0 3 j 


FCO? N P 1 95 

— - — — + Die 

p p S 0 at 

(B15) 


B8 Alkalinity 


dTA 

d? 



- (ta* - ta) 


. N P 1 95 
TA) = — + - — TA 
p S 0 at 


CaC0 3 


CaC0 3 ■ 


PcaC0 3 ] 


(B16) 


B9 Calcite production, calcification/organic carbon 
production ratio (RcaC 03 ) and calcite dissolution 
rate (Xc a co 3 ) 


T’CaCOs = ^CaC0 3 
[(/A3,phyto pi 0.5£ ZO o,NO 3 ^o) 3Tphy,det( P P.o)] 


Coc 


(B17a) 


^CaC0 3 = ^CaCOp F 3 , p h y (£) min 

Niim, — 1 max \ 0.0001 , 

AT 1 +PO 4 


2 + T 


1.*' 

2 


^CaC0 3 — ^CaC0 3 


ACO 3 

*CaC03 + ACO 3 
ACO 3 = max(0, CO 2 .^ — CO 2- ) 

B10 Net community production 


(B17b) 


HI 


I [/^l.phy^l “1“ M2,phy^2 (^p^phy,det “1“ ^phy,N0 3 ) 

/ Red 

{(P\ - P\,o) + (Pi - P 2 ,o)}~a z M zooM Z 2 (B18) 

^zoo,no 3 (Z Zq) rem det^J k TC DOC+ ^ M 3 ,phy P 


Bll Oxygen 


dO 2 FO2 Np / O2 \ 

d t Z P P V N J Red 

B12 Nutrient limitation 


(B19) 


Ni im = NH 4lim +N0 3llm 


(B20) 


N°3 (I-NH 4 ) 

N03lim = 

(An03 + NO 3 ) (*NH4 + NH 4 ) 

(B21) 

nh 4 

NH4i im = 

(A NH4 +NH 4 ) 

(B22) 

NH 4 ij m 

NH 4 ii m + NC>3ii m 

(B23) 

NC>3ii m 

NH 4 li m + N03i; m 

(B24) 


B13 Irradiance model 

The total (infrared plus visible) solar radiation is obtained us- 
ing the Frouin model (Frouin et al., 1989). This model pro- 
vides the total radiation and the photosynthetically available 
radiation (PAR). The infrared (/ir) component is obtained by 
subtracting the PAR component from the total solar radiation. 
Using a spectral model for PAR (Gregg and Carder, 1990), 
the spectral PAR c o in p o 11 c n t / 1 > A R ( }, jean be determined. The 
infrared component (for mixed layer model only) and the 
PAR component (mixed layer and biogeochemical model) of 
the penetrating irradiance are obtained from 

/ir(z) = /ir(z - Az)exp[— « irAs], (B25) 


fpAR (Pz) = /(A.,z- As)exp[-(a u ,(A) + fl ph (A))Az] (B26) 

where antO^m' 1 ) is the attenuation coefficient for in- 
frared radiation and a w (X) and fl p h(l) are the wavelength- 
dependent light attenuation coefficients for water and phyto- 
plankton, respectively. The water and chlorophyll-dependent 
attenuation coefficients from Morel (1988) were used in the 
model for this study. The dissolved matter attenuation coef- 
ficients, adm(X), are calculated by applying the IOP (inherent 
optical properties) model of Garver and Siegel (1997), which 
uses water leaving radiances from 6 SeaWiFS bands as in- 
put (level 3 binned monthly composites). The IOP model 
calculates the attenuation coefficient due to dissolved matter 
for the 443 nm wavelength. The attenuation coefficients for 
other wavelengths are obtained from 


fldm(A-) = fldm(443)exp[S(A.-443)] (B27) 


where the exponential decay constant, S, is chosen to be 
0.02061. A correction is applied to the downward irradiance 
pathway to account for seawater light refraction following 
Snell’s law. After some algebraic manipulations, the correc- 
tion is applied toAzas follows: 


A z' = 


Az 
cos /I 


(B28) 


P = a sin 



(B29) 


where a z is the solar zenith angle, and n s is the seawater re- 
fraction coefficient, which is expressed as a function of salin- 
ity and temperature adapted from Table 3.12 of Neumann and 
Pierson (1966) as 

n s = 10“ 6 (285.77 - 15.65T + 197.675) + 1.333338 (B30) 
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B13.1 Ammonium nitrification 


B18 Model forcing and relaxation approach 


An = A“ aX I 1 - 


r=24 A=470 


D-D n 


D — I) mm ~ Kd 


°-I / 


I ( X,t)as (k)dXdt 

t= 0 A=300 

B14 Phytoplankton growth 


(B31) 


M/,phy(N03,P0 4 , Si0 2 , E) = F it p hy (£) 
P0 4 Si0 2 


Nln 


£po 4 + PO4 k$io 2 + Si0 2 J 


(B32) 


M”phv«i 7 

F /iPhy (/) = ' P ^ = (B33) 

v (<&) 2+( “'- /)2 

Note: Si0 2 included for diatoms only. ri=14 is N:P. 

B15 Zooplankton grazing 

P 2 

P*oo(Pi ) = /Co X 7 (B34) 

£phy + Pf 

B16 Primary production 

PP = [(/^l,phyPl + M2,phy^°2) ( ~ ) + M3,phy^3 ( ~ ) ]^MW 

V^/ Red \ w / Coc 

(§L = 6 - 625 (§L = 9 - 4 c ™= 12 (B35) 

( Si\ _ 15 

\ ^ / Red 16 

The C:N ratio for coccolithophores (Coe) of 9.4 is the av- 
erage from the reported range of 5.81 to 13.05 in Fernandez 
et al. (1993). 

B17 Chlorophyll and Chl:N ratio 

The model chlorophyll is calculated following the photoad- 
aptation scheme for Chl:N ratio of Doney et al. (1996): 

Chla=Chl:N^P/ (B36) 


Chl:N=Chl:N max - (Chl:N max - Chl:N min ) ^ / PA r < /* 
Chl:N=Chl:N min / PAR > /* 

(B37) 

The subscripts phy, zoo, and det refer to phytoplankton, 
zooplankton, and detritus, respectively. Table B 1 defines the 
model state variables and Table B2 provides the definition 
of the ECO1D-2.0 parameters and values used. Equivalent 
parameters for ECO ID- 1.0, configured for and applied at the 
BATS site, are shown in parenthesis for comparison. 


The terms <5(z)FC0 2 /p and <5(z)FO 2 /p, in (B15) and (B19), 
respectively, represent the C0 2 and 0 2 sea-air fluxes at the 
surface. The Kroenecker delta (5[z=0]=l ; <5[= >0]=0) is used 
to denote that carbon dioxide and oxygen fluxes (FC0 2 and 
F0 2 , respectively) are only applied at the sea-air interface. 
The following formulations for the C0 2 and 0 2 gas transfer 
were applied in the form of flux boundary conditions (FC0 2 
and F0 2 in mmol m -2 yr -1 ) at the sea-air interface: 

FC0 2 = K 0 a A /2 C0 2 F0 2 = K a [O? - 0 2 ] (B38) 

where, K a is the gas transfer velocity, in m d -1 , which 
is a function of water temperature and wind speed (Wan- 
ninkhof, 1992), a is the C0 2 solubility in seawater (in mmol 
m -3 /x atm), which is a function of temperature and salinity 
(Weiss, 1974), ApC0 2 (in /x atm) is the difference between 
sea and air pC0 2 , and O* is the oxygen saturation concen- 
tration (in mmol m -3 ) in seawater, which is a function of 
temperature and atmospheric pressure (Weiss, 1970). 

We adopt the following relationship between gas transfer 
and wind speed ( W ) (Wanninkhof, 1992) using the NCEP 3 
hourly winds: 

K a — 0.31VT 2 (Sc/660) -1 / 2 (B39) 

where Sc is the Schmidt number of C0 2 or 0 2 (Wanninkhof, 
1992). 

To account for horizontal advective processes of heat and 
salt within deeper layers of the 1-D mixed layer model, tem- 
perature and salinity are assimilated from the 3-D model us- 
ing a straightforward approach. The approach consists of re- 
laxing the temperature and salinity profiles calculated by the 
1-D mixed layer model to the pre-calculated values provided 
by the 3-D model at all depths. The assimilation of T and 
S is done using a Newtonian relaxation (nudging) method 
(Bauer and Wulfmeyer, 2009) with a relaxation time scale 
(r) of 10 days for T and 30 days for S. For model proper- 
ties where the relaxation time scales for deep water are short 
(10 days or less), the term “nudging” is not applicable, as the 
properties are set to observed values. 

A similar relaxation approach is used for nitrate, phos- 
phate, and silicate, except that the nutrient values originate 
from '/’-dependent equations obtained from T, NO3, PO4, 
and Si0 2 climatologic monthly profiles (0-500 m) from the 
World Ocean Atlas 2005 at the model site. The total num- 
ber of data points is 168 (/V=14 depthsxl2 months). The 
equations are 

NO^ = -2.253(±0.169)7’ + 29.92(±1.21) 
r 2 = 0.814 RMSE= 1.258 
PO4 = — 0. 1333 (±0. 0095) T + 1.884(±0.068) 
r 2 = 0.836 RMSE = 0.0707 
SiO, = —1.479(3=0.135)7’ + 17.37(3=0.98) 
r 2 = 0.735 RMSE =1.012 (B40) 
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DIC* = nC T x 5/35 

nC T = 1962.6- 14.6204(T- 20) -0.1371(7’ -20) 2 (B41) 

r 2 = 0.762, RMSE = 10.15, N = 148 

TA* =231 1.6 + 46.4153(5' -35) + 56.4425(5 - 35) 2 (B42) 
-0.0456(7 - 20) - 0.0387(7 - 20) 2 
r 2 = 0.448, RMSE = 4.15, N = 123 
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